Classical Nucleation Theory of the One-Component Plasma 
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We investigate the crystallization rate of a one-component plasma (OCP) in the context of classical 
nucleation theory. From our derivation of the free energy of an arbitrary distribution of solid clusters 
embedded in a liquid phase, we derive the steady-state nucleation rate of an OCP as a function of 
the Coulomb coupling parameter F. Our result for the rate is in accord with recent molecular 
dynamics simulations, but it is greater than that of previous analytical estimates by many orders of 
magnitude. Further molecular dynamics simulations of the nucleation rate of a supercooled liquid 
OCP for several values of F would clarify the physics of this process. 



I. INTRODUCTION 

The one-component plasma (OCP) consists of a sin- 
gle species of N charged ions of mass to in a volume V 
in a uniform neutralizing electron background. This ide- 
alized system is a model for astrophysical settings such 
as neutron star crusts and white dwarf interiors. The 
dimensionless coupling parameter 

characterizes the state of the OCP, where Ze is the ion's 
charge, a — (3/47mi)^/'^ is the mean distance between 
ions at number density ni — N/V, and T is the tem- 
perature. We consider the ions to be classical with a 
deBroglie wavelength A = (h'^ /27rmkBTy^^ ^ a. For 
r ^ 1, the Coulomb coupling is weak and the ions be- 
have like an ideal Maxwell-Boltzmann gas. For F ^ 1, 
the Coulomb coupling is large and the system is in the 
liquid phase until F > F,„ « 175 [lj,l3i when the ions un- 
dergo a first-order phase transition into a periodic lattice 
with a body-centered cubic configuration. 

The Helmholtz free energy of the OCP in the liquid 
phase is written as F = F"* -I- F°^, where the ideal part 
F"* — Nk^T [ln(7iiA'^) — l] , and interactions are incor- 
porated into the 'excess' part F°^ = —kBTlog{ZN/V^), 
where Z^v is the configuration integral [3|. In the non- 
interacting limit, Zn -^ V^ , and F''^ — > 0. Simula- 
tions of the liquid OCP have generated accurate fitting 
formulas for the reduced excess free energy /""(F) = 
F^^/Nk^T [e.g., Jj|, and calculations of the solid state 
free energy [e.g.,l5| allow for the determination of F^. 

Crystallization of an OCP occurs via nucleation, the 
process by which small crystals form via localized fluc- 
tuations in the liquid and grow via the accretion of sur- 
rounding ions [see e.g.,|6[. Empirically, most pure liquids 
can be supercooled (in our case F > Fm) without solidi- 
fying, which implies that a kinetic barrier to nucleation 
exists [e.g., lD]. Monte Carlo [a, l9| and molecular dy- 
namics jlO| simulations demonstrate that such a kinetic 
barrier exists in an OCP as well. Localized fluctuations 



of sufficient amplitude can generate crystals large enough 
to overcome the kinetic barrier, and hence they are ther- 
modynamically stable. Such crystals subsequently grow 
and facilitate the phase transition. 

The rate at which astrophysical solids grow can af- 
fect the degree of chemical phase separation during the 
solidification of both white dwarf interiors [ll|, and refer- 
ences therein] and neutron star crusts |1J| , the existence 
of an amorphous glassy state in white dwarf interiors 
131 Il4 |. and the prevalence of defects and impurities in 
neutron star crusts [15|, [l^]. Furthermore, the time at 
which the latent heat of crystallization is released within 
a white dwarf affects its observed cooling rate [e.g., llTI . 
and references therein]. Unfortunately, the nucleation 
rate of an OCP is not well understood. Several authors 
[l4 llSl . llSl . Iig ] assumed that classical nucleation the- 
ory [6^ applied to the OCP, but the accuracy of their 
assumptions was untested. Daligault [lOJ was the first 
to directly study nucleation of an OCP using molecular 
dynamics simulations. Although his results were qualita- 
tively consistent with classical nucleation theory, it was 
unclear whether or not the results agreed quantitatively. 
In this paper, we reexamine the nucleation of an OCP in 
the context of classical nucleation theory, but with a dif- 
ferent approach to the statistical physics of an ensemble 
of solid clusters embedded in a liquid phase. We find a 
steady-state nucleation rate as a function of F in quanti- 
tative accord with those of flOj. 

We begin in Sec.|TT]by reviewing what is known about 
the excess free energies for multi-component plasmas and 
then deriving the minimum reversible work of cluster for- 
mation and the equilibrium distribution of cluster sizes. 
We use these results to deduce the nucleation rate in 
Sec. mil and we compare this rate, along with that used 
in previous work, to the results of [lO| in Sec. IIVI We 
close in Sec. |V] by discussing our results and suggesting 
future numerical experiments and analytic improvements 
to this work. 



II. MINIMUM REVERSIBLE WORK OF 
CLUSTER FORMATION 

We consider a classical OCP consisting of Nj ions of 
mass m in a uniform neutralizing background of electrons 
at Ue = Zn\. Motivated by the typical astrophysical con- 
text, we also assume that degenerate electrons supply all 
of the pressure so that rie is constant during any phase 
change and set by charge neutrality. However, before 
discussing how we treat solid clusters embedded in a liq- 
uid OCP, it is important to clarify what is known about 
mixtures of ions of different charges. 

Hansen tt al. |2(l| first showed that the excess 
Helmholtz free energy in a mixture of ions of charges 
Zi and number Ni in a volume V is simply given by the 
sums of the excess free energies of each separate species 
at the same electron density, 



F'^'^ ^ kBTy N,mr,) 



(2) 



Z,- Fe and Te — e^ / aeksT , where Oe 



Here Ti = 

(3/47rno)^^'^ is the mean distance between electrons, and 
all Fi > 1. Recent simulations [^ have shown that the 
deviations from this simple linear mixing rule are less 
than 0.05%. Hence, to an excellent approximation, the 
free energy of a mixture of liquid ions is simply 
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spherical, and that all clusters of this size conform to this 
configuration. 

The total Helmholtz free energy of the system is 
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where fik^T and fsk^T are the free energies of one ion in 
the pure liquid and pure solid state, respectively, aksT 
is proportional to the surface free energy of an ion at 
the boundary of a cluster, Ni and Zi = iZ are the num- 
ber and charge, respectively, of clusters (including liquid 
monomers) that consist of i ions. 



N, 



total 






(6) 



is the total number of clusters (including monomers), and 

N,. 



{z) = Y.z, 



N, 



total 



is the average charge of the clusters. Equation 
subject to the constraint 



^i = E*^- 



(7) 



IS 



(8) 



Later authors chose to write this equation differently by 
using the fits, /;, of the total Helmholtz free energies of 
each pure state and then mixing them. This results in 
the more often seen equation 



El2l^ = y N,fi{T,zl") + y TV, In (^^ 
knT ^ ■' ^ ' "^ \{Z)Ni 



total 



where (Z) is the average charge. The last term is minus 



2ll and the Appendix] , 



the entropy of mixing [see e.g.,|2C 
although it originates simply from charge conservation at 
a fixed Ue- This term will appear in our derivation of the 
minimum reversible work of cluster formation. 

This knowledge of mixing at fixed rig motivates our 
re-derivation of the minimum reversible work of cluster 
formation, and thereby the equilibrium cluster distribu- 
tion, from the total Helmholtz free energy i^totai of an 
arbitrary distribution of clusters. Following Frenkel [23 
(see also [23|,|24|), we consider clusters that contain dif- 
ferent numbers of ions to be distinct species, and we con- 
sider all clusters that contain the same number of ions 
to be indistinguishable. The latter is an approximation 
because such clusters could have different morphologies. 
Thus we assume that, for each cluster size, there exists 
a unique configuration of ions for which the total sur- 
face free energy is a minimum, that this configuration is 



From equation ^, the chemical potentials of 
monomers and clusters are 
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respectively. Equation ©, the chemical potential of a 
liquid monomer, is identical to that of [25]. According 
to classical nucleation theory, clusters grow or shrink by 
the addition or subtraction of single monomers. Thus in 
chemical equilibrium [li = /ii_i + /ii, and so 



[li ^ ifii 



fll) 



by iteration. From equations (O, (fTU|) . and pTj) . and not- 
ing that (.^)A^totai = ZNi, the equilibrium distribution of 
clusters is 



iN, 



^=(^^yexp(-[z(/.-/,-l) + zV3^ + l]), 



(12) 

which is in the form of the law of mass action. The 
monomer concentration, and hence the concentration of 
all clusters, is found by invoking the constraint of equa- 
tion ([5]). Alternatively, equation (|12p can be written in 



terms of the Boltzmann distribution P{i) — exp(— A/i), 
where P{i) = iNi/Ni is the probabihty that an ion picked 
at random is a member of a cluster of size i and 



A/, - 



-hi(7Vi/7Vi), 



1=1, 
a + 1, i>2 

(13) 
is the minimum reversible work divided by k^T of form- 
ing a cluster of size i. Demanding that P{i) is normalized 
then yields Ni/Ni. 



III. STEADY-STATE NUCLEATION RATE 

We follow the treatment of classical nucleation theory 
by Kelton ^6] and derive the steady-state cluster nucle- 
ation rate. We consider an OCP at fixed F that is entirely 
in the liquid state at time i = and calculate the rate of 
cluster formation for t > 0, when A^i w iVi, and so 



Af, = i{f,-fi-l) + i 



2/3^ 



1, 



(14) 



from equation (fT3|l . We use the fit of [1] for fi{T), where 
we take the coefficients derived from the simulation data 
of [26|], and we use the fit of [S] for /s(r). For sufficiently 
large F, (/s — /; — 1) < 0. The fit of fg given by equations 
(10) and (15) of Q is vahd only for 170 < F < 2000, and 
the above inequality is easily satisfied for the entire range 
of F in which the fits are applicable. The minimum F 
above which (/g — /; — 1) < is unknown. Since a > 0, 
there exists a critical size 



2a/3 



l + fi-fs 



(15) 



such that the reversible work of cluster formation, A/^, 
is a maximum. Cluster growth is energetically favorable 
for clusters of size i > i* and unfavorable for clusters of 
size i < i*. Equation (|14p can now be written in the 
convenient form 
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-21- 



+ 1. (16) 



In classical nucleation theory, crystallization occurs via 
the formation and subsequent growth of critical clusters 
of size i*. Clusters of size i < i* are formed via equi- 
librium fiuctuations and are presumed to be transient. 
Equation (fT2|) gives their distribution. Clusters of size 
I > «* are presumed to be stable. Such clusters only grow 
with time. The formation of stable clusters is therefore a 
two-step process: a stable cluster forms when (1) equilib- 
rium fluctuations in the liquid phase generate a transient 
cluster of size i*, and (2) that cluster accretes an addi- 
tional monomer. Therefore, the steady-state nucleation 
rate J, defined to be the number of solid clusters formed 
per second, is roughly equal to the number of critical 



clusters of size i* times the rate at which a monomer at- 
taches to a critical cluster. Kelton [6| performs a more 
thorough derivation of the nucleation rate and finds 



J(7Vi,F) = 



24:Di 



.2/3 



-N^ 






1/2 



exp(-A/,J, (17) 



where D is the diffusion coefficient and (A/i^/37ri^)^/^ 
is the Zeldovich factor. We use equation (fT7|) for the 
steady-state nucleation rate of an OCP. For D, we use 
the classical OCP diffusion coefficient [27(| 



D = 2.95u;„a^r 



^1.34 



(18) 



where Wp = [47mi(Ze)^/rn]^'^ is the ion plasma fre- 
quency. Recent calculations have verified the accuracy of 
this expression for F < Fn 



28|, l29| as well as for F > Fn 



. For the normalized surface free energy a, we follow 
and [ISl and set 



a = 4:71 



2/3 



4^ j ^^ kBT 



^^^"^ = {36ny/'(3\ ( -^ ) , (19) 



F 



where /3 « 0.5 is the ratio of the interfacial energy and 
heat of fusion empirically derived from liquid metal ex- 
periments [e.g., |3l|, I34I, and Afce^^m is the enthalpy of 
melting per ion. For this work, we set A = 0.77 [5|, |33|. 

IV. COMPARISON TO PREVIOUS WORK 

Several authors have investigated the homogeneous nu- 
cleation rate of an OCP [14i, [la, |l8|, ll9|. However, in 
deriving the minimum reversible work of cluster forma- 
tion A/i, they assumed that the entropy of mixing is 
simply that of an ideal gas, i.e., the term hi{Zi/{Z)) 
in equation ([5]) was omitted. Consequently, they find 
A/i = i(/s - fi) + i'^^^a, and equations ([TS]) and (fTB]) 
become 






A/, 



2a/3 
/; - /s 

(W3)^ 
2{fi - fsf 



2/3 



(20) 

(21) 



We demonstrate below that previous calculations result 
in a nucleation rate that is orders of magnitude lower that 
the rates deduced from simulations, whereas our model 
predicts a rate closer to that found numerically. The 
steady-state nucleation rate per ion as a function of F for 
both models is shown in Figure 1. 

Presently, crystallization of an OCP is studied only via 
numerical simulations. Although there are many calcu- 
lations of the free energy difference of the pure liquid 
and pure solid states as a function of F, few have in- 
vestigated the onset of crystallization. Ogata [8] and De- 
Witt et al. [9| studied nucleation in an OCP using Monte 
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FIG. 1: Steady-state nucleation rate per ion as a function of 
r. "This Work" is the rate calculated from our model, and 
"Previous Work" is the rate calculated using the models of 
previous authors. Our calculations predict that the nucleation 
rate of a one-component plasma is > 10''' times higher than 
previously assumed for all values of V . 



Carlo simulations. Unfortunately, time is not a variable 
in Monte Carlo simulations, so that the rate of nucleation 
is difficult to calculate [although see [34 1 . 

Recently, Daligault [lO| performed molecular dynam- 
ics simulations to investigate the kinetics of nucleation in 
a supercooled OCP. Using N\ = 4394 ions, he calculated 
the total number of solid nuclei as a function of time 
for r = 400 and 300. In the V = 400 calculation (see 
Figure 4 of [101), following a short period of transient 
nucleation, there is a significant period of time during 
which steady-state nucleation takes place, as evidenced 
by the nearly constant slope in the plot of the number 
of clusters as a function of time. In the F — 300 calcula- 
tion (see Figure 8 of [l0|), there is no obvious period of 
time during which steady-state nucleation occurs. To be 
consistent with our estimate from the previous calcula- 
tion, we estimate Jmd, the steady-state nucleation rate 
deduced from molecular dynamics simulations, from the 
period of time during which the nucleation rate is nearly 
constant, which occurs just before nucleation saturates. 
In Table I, we list the inferred values of Jmd as well as 
the nucleation rates calculated using both equations (jl51 
[T5)) of our model and equations (|20ll2ip from previous 
models. Our model predicts nucleation rates that are in 
reasonable agreement with the results of Daligault [lO|, 
whereas the predictions of previous models are far too 
low. 



V. DISCUSSION 

Through an improved treatment of the statistical 
physics of the distribution of solid clusters in the liq- 



TABLE 1: Steady-state nucleation rates and critical cluster 
sizes for Ni = 4394: "Jmd" is the inferred nucleation rate 
from the simulations of [l(]( |. "(this work)" denotes a quan- 
tity calculated from our model, and "(previous)" denotes a 
quantity calculated using the model of previous authors. 



r Jmd/ia'p 


J/a;p J/tJp i, i» 




(this work) (previous) (this work) (previous) 


400 ^ 0.2 


0.4 1 X lO"'' 4 51 


300 ~ 0.1 


2 2 X 10"* 3 91 



uid phase, we have developed a model of steady-state 
nucleation in an OCP, and the results of our model are 
in accord with those of the time-dependent molecular dy- 
namics simulations of 'lCf|. Our results suggest that crys- 
tallization in an OCP occurs at a rate that is more than 
five orders of magnitude higher than previously assumed. 
Nucleation in an OCP differs from nucleation in many 
other media because of two unique features of an OCP: 

(1) The degenerate electron background, which supplies 
the pressure, is unaffected during the phase transition. 
Consequently, both pressure and volume remain con- 
stant. Although the complete non-interaction of the elec- 
trons is an approximation of the OCP model, the role 
electrons play in more realistic models is negligible [l|]. 

(2) The binding of ions of like charge in a crystal is weak. 
This is reflected in the difference of the liquid and solid 
free energies. Indeed, this difference does not exceed k^T 
(i.e. Ji-U does not exceed 1) until F/Fm = Tm/T > 2.8. 
The combination of these two features has interesting im- 
plications for the crystallization of an OCP. In particular, 
nucleation in an OCP is driven to a significant degree sim- 
ply by the increase in entropy resulting from cluster for- 
mation: the change in the total energy of a supercooled 
OCP when a cluster forms is small relative to the in- 
crease in the entropy. In fact, equation p4p implies that 
a population of solid clusters will exist at a given time 
for a range of F < Fm because the entropy of an ensem- 
ble of clusters embedded in a liquid OCP is higher than 
that of the pure liquid phase. This may be related to the 
caging effect observed in molecular dynamics simulations 
of an OCP in the liquid phase 28J, |35|. Furthermore, we 
note that a reanalysis of the entropy of mixing term and 
thereby the classical nucleation rate may be needed for 
other supercooled liquids. 

There are a few important issues that we have not ad- 
dressed. The surface free energy of a solid cluster in a 
Coulomb liquid is not known. The expression for the 
surface free energy contribution used in this work is only 
our best estimate, but we cannot be confident that the 
numerical coefficient in equation (|19p is accurate. This 
issue requires further investigation via molecular dynam- 
ics simulations. Some of the assumptions we made in 
Sec. In] may not be valid. We assume that solid clusters 
grow or shrink by the addition or subtraction of liquid 



monomers, but the simulations of Daligault [10| imply 
both that two solid clusters can fuse together and that a 
solid cluster can fission into two or more clusters. Fur- 
thermore, we assume that all clusters of a given size have 
the same morphology and are thus indistinguishable, but 
the simulations of Daligault [lO| show that two clusters 
of the same size can in fact have different morphologies. 
We currently do not know what effects, if any, a violation 
of either of these assumptions would have on our results. 
Further molecular dynamics simulations are impera- 
tive to understand nucleation in an OCP. It is evident 
from equations (|17m9p that the steady-state nucleation 
rate may be written as 



J = NiLUpF{r), 



(22) 



where F(T) is a function only of F. Numerically con- 
firming that the nucleation rate scales linearly with Ni 
and LUp should be simple. One could then deduce -F(F), 
and hence J, by measuring the nucleation rate for sev- 
eral different values of F. Calculations for F^ ^ F < 210 
would be particularly useful. Figure 1 shows that there is 
an enormous difference between our results and those of 
previous authors: our model predicts that some clusters 
should form rather quickly, whereas other models predict 
that clusters essentially never form. 
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APPENDIX: ENTROPY OF MIXING IN A 
ONE-COMPONENT PLASMA 



We follow the treatment of [36| and derive the ideal 
entropy of mixing of clusters in an OCP. The entropy of 
an ideal system of Ni indistinguishable particles of mass 
mi at temperature T and within a volume V is 
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The ideal entropy of mixing is thus 
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Admixing = 2_^ [<5'i(Vtotal) — SiiVi)] , 
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is the fractional volume occupied by the Ni clusters of 
charge Zi and 
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is the total volume of the system. The Vi's are set by 
demanding both local charge neutrality and that ric is 
constant throughout the medium. It follows that 
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and hence 
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